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Introduction: Numerical modeling of impact cra- 
tering has reached a high degree of sophistication; 
however, the treatment of porous materials still poses a 
large problem in hydrocode calculations. Porosity 
plays only a minor role in the formation of large cra- 
ters on most planetary objects, but impacts on comets 
are believed to be highly affected by the presence of 
porosity, which may be as much as 80%. The upcom- 
ing Deep Impact Mission (launched January 2005) will 
provide more detailed data about the compostion of a 
comet (Tempel 1) by shooting a -370 kg projectile 
onto the surface of its nucleus. The numerical simula- 
tions of such impact events requires an appropriate 
model for how pore space in the comet is crushed out 
during the violent initial stage of the impact event. 

Most hydrocodes compute the pressure explicitly 
using an “equation of state” (EOS) for each material, 
which relates changes in density and internal energy to 
changes in pressure. The added complication intro- 
duced by porosity is that changes in a material’s den- 
sity are due to both the closing of pore space (< compac- 
tion ) and compression of the matrix. The amount of 
resistance to volume change and the amount of irre- 
versible work done during these two processes is very 
different; it is far easier to compact a porous material 
sample than to compress a non-porous sample of the 
same material. 

As an alternative to existing porosity models, like 
the P-a model [1], we present a novel approach for 
dealing with the compaction of porosity in hydrocode 
calcualtions. 

Porosity Model: Any initial porosity in the target 
results in an increase in irreversible work during the 
early stages of the impact due to the collapse of pore 
space and the corresponding increase in density at rela- 
tively low shock pressures. Thus higher temperatures 
occur in the material after the shock wave has passed 
and, consequently, the Hugoniot curve of pressure 
versus specific volume is much less steep for the po- 
rous sample (see Fig 2). As the pore space is reduced, 
the slope of the pressure-specific volume curve steep- 
ens; the slope becomes identical to that of the non- 
porous material when all pores are closed. Once all 
pores are removed, however, the Hugoniot curves of 
the porous and non-porous sample of the same mate- 
rial do not lie on top of one another. This is because a 
greater amount of irreversible work has been per- 
formed on the porous sample, raising its internal en- 


ergy to a higher level and hence raising the pressure at 
an equivalent specific volume. 

Previous models, like the often used P-a model [1], 
account for the collapse of pores by assuming a pres- 
sure-dependent compaction function of the porosity (/). 
Instead of porosity, a distention parameter a>l is often 
used, with a= 1/(1 -(/)), a=V/Vs (V specific volume, V s 
specific volume of the matrix). The pressure P is then 
calculated as follows [2]: 


P = (E, p,a) = — P s (< ap, E) = —P s (p s , E) 
a a 


(i) 


where P s is the pressure in the solid matrix, p s the den- 
sity of the matrix, p bulk density, and E is the internal 
energy. This equation causes some problems when 
implemented into a hydrocode: For a given time step t , 
a t + ; must be known to derive P t+1 , but a t+1 =f(P t+1 ). A 
common solution is to use small sub-cycles to iterate 
the new P t+1 value. This method requires extra compu- 
tation time and is numerically not very stable. 

As an alternative, we propose to use a volumetric 
strain dependency for the compaction function 
a=f(£ v ), instead of pressure. s v is more directly related 
to the crushing of pore space than pressure and a 
maximum compaction rate can be readily derived from 
the definition of volumetric strain: 
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If one assumes that no compression of the matrix oc- 
curs during the crushing of pores (F,=constant) eq. 2 
yields: 
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In reality the crushing-out and closure of pores 
does not begin until a critical volumetric strain is 
reached s e . This defines a kind of plastic yield thresh- 
old below which volume changes are elastic (reversi- 
ble) and above which compaction (permanent pore 
collapse) occurs (see Fig. 1). Moreover, pore collapse 
may not occur at the maximum rate and some com- 
pression of the matrix accompanies the compaction. 
Hence, we introduce a parameter k<1 to control the 
rate of pore closure in the compaction regime. The 
exponential compaction (eq. 2) works well for the 
early stage of compaction, but does not give a smooth 
transition as a approaches 1 (complete compaction). 
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As a result we define a secondary compaction regime 
for volumetric strains less than the transition volumet- 
ric strain s x (£v <8 x)- We use a power-law expression 
for the compaction rate (as is frequently used in im- 
plementations of the P-a model): 


a = f(e v ) = a 0 e K{£v ~ £e) 
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Eq. 4 applies for volumetric strains greater than s c , 
the volumetric strain where all pore space is crushed 
out; further volumetric strains are purely compressive. 
Note, that s c can be derived for a given s x to achieve a 
continuous transition at a=f(£ x )='F(£ x ); thus, only s & 
s x and /rare free parameters. 



Fig 1: Compaction model showing the volumetric strain 
dependence of the distention a. Note that in this figure posi- 
tive strain is compressive whereas in all other equations in 
this paper s v >0 is tensile. 

Verification of the model: The present approach 
was incorporated in the SALE-3 MAT hydrocode [3], 
which is based on the SALE hydrocode [4] and is ba- 
sically very similar to SALEB by Ivanonv et al. [5]. To 
verify the model we first (i) reproduced data from 
shock experiments on porous Aluminum for both Eul- 
erian and Lagrangian mode. Fig 2 shows modeled 
Hugoniot curves in comparison to experimental data 
[6] with porosities (p=0 ( a=l ) and <fi=0.5 ( a=2 ). A 
good fit was obtained by adjusting the free parameters 
for the porosity model. We used a Tillotson equation 
of state to calculate pressure from density and internal 
energy (see eq. 1). (ii) Additionally, we modeled a 
simple one-dimensional impact experiment (reported 
by Butcher et al [7]) in which a porous Aluminum 
sample ( a=1.3 ) was impacted by a thick Aluminum 
plate, with a velocity of 126 m/s, and buffered by an- 
other Aluminum plate. A velocity interferometer 
(VISAR) was used to measure the velocity at the free 
surface of the buffer plate as a function of time. Fig 3 


shows the results of the numerical model in compari- 
son to the experimental data. 
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Fig 2: Calculated Hugoniot curves for Aluminum with a=l 
and a=2 in comparison to experimental results [2]. 



experiment as described in the text. 

Conclusion: We present a simple and intuitive po- 
rosity model for use in hydrocodes that reproduces 
experimental data well. The model requires only three 
material- specific input parameters and each has a 
physical basis. Thus, although these values should be 
determined for a given material by comparing model 
results with experimental data, it is possible to make 
sensible estimates for the parameters where no data for 
a porous material exists. For future work we intend to 
transform the two-regime compaction definition (eq. 4) 
into one equation that fulfills the given boundary con- 
ditions. 
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